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ABSTRACT 

We study the impact of baryonic physics on cosmological parameter estimation with weak lensing 
surveys. We run a set of cosmological hydrodynamics simulations with different galaxy formation 
models. We then perform ray-tracing simulations through the total matter density field to generate 
100 independent convergence maps of 25 deg 2 field-of-view, and use them to examine the ability 
of the following three lensing statistics as cosmological probes; power spectrum, peak counts, and 
Minkowski functionals. For the upcoming wide-held observations such as Subaru Hyper Suprime- 
Cam (HSC) survey with a sky coverage of 1400 deg 2 , these three statistics provide tight constraints 
on the matter density, density fluctuation amplitude, and dark energy equation of state, but parameter 
bias is induced by the baryonic processes such as gas cooling and stellar feedback. When we use power 
spectrum, peak counts, and Minkowski functionals, the magnitude of relative bias in the dark energy 
equation of state parameter w is at a level of, respectively, Sw ~ 0.017, 0.061, and 0.0011. For HSC 
survey, these values are smaller than the statistical errors estimated from Fisher analysis. The bias 
can be significant when the statistical errors become small in future observations with a much larger 
survey area. We find the bias is induced in different directions in the parameter space depending on the 
statistics employed. While the two-point statistic, i.e. power spectrum, yields robust results against 
baryonic effects, the overall constraining power is weak compared with peak counts and Minkowski 
functionals. On the other hand, using one of peak counts or Minkowski functionals, or combined 
analysis with multiple statistics, results in biased parameter estimate. The bias can be as large as 
lfj for HSC survey, and will be more significant for upcoming wider area surveys. We suggest to use 
an optimized combination so that the baryonic effects on parameter estimation are mitigated. Such 
‘calibrated’ combination can place stringent and robust constraints on cosmological parameters. 

Subject headings: gravitational lensing: weak — cosmological parameters 

cosmology: theory - large-scale structure of the universe 


1. INTRODUCTION 


An array of recent observations of the large-scale struc¬ 
ture of the universe such as cos mic microwave back¬ 
ground (CMB) anisotropies (e.g. , |Hinshaw et al.||2013 
Planck Collaboration et al. 2014) and galax y clustering 


(e.g., Reid et al. 2011) "Beutier et al. 2014) established 
the standard cosmological model called AGDM model. 
In ACDM model, the energy content of the present-day 
universe is dominated by two mysterious components: 
dark energy and dark matter. Dark energy realizes the 
cosmic acceleration at present and dark matter plays an 
important role of formation of rich structure in the uni¬ 
verse. However, we have not understood yet the nature 
of dark energy and the physical properties of dark mat¬ 
ter. In order to reveal the mysterious dark components 
in the universe, several observational programs are pro¬ 
posed and still under investigation. Such observational 
programs include Subaru Hyper Suprime-Cam (HSCQ 
the Dark Energy Survey (DES^ and the Large Synop¬ 
tic Survey Telescope (LSST^ Space missions such as 


' E-mail: ken.osato@utap.phys.s.u-tokyo.ac.jp 

1 http://www.naoj.org/Projects/HSC/index.html 

2 http://www.darkenergysurvey.org/ 

3 http://www.lsst.org/lsst/ 


EuclicF] and WFIRSfj^] are also promising. Gravita¬ 
tional lensing is expected to be the main subject of these 
future surveys that are aimed at studying the large-scale 
structure of the universe at present and in the past. 

Weak gravitational lensing (WL) by large-scale struc¬ 
ture in the universe is the promising probe into proper- 
ties of dark matter and dark energy (for a review , see 
Bartelmann fe; Schneider||2001[ IMunshi et al.|[2008l |Kil-| 
binger 2014[). WL causes small distortion of image of 


distant source galaxies, called cosmic shear, which reflect 
directly the intervening matter distribution along a line 
of sight. Recent cosmic shear observations have proved 
WL measurement to be a powerful tool for studying dark 
matter distribution in the universe, from which one can 
extract in formation on the ba s ic cosmological param e- 
ters (e.g., |Massey et al.| (|2007|); |Kilbinger et al.| (|2013|)). 
Forthcoming weak-iensing surveys are aimed at measur¬ 
ing cosmic shear over a wide area of more than 1000 


deg 2 . These observations will address important ques¬ 
tions of dark matter and dark energy at unprecedented 
precision. 

Unfortunately, major statistical methods to make the 


4 http://www.euclid-ec.org/ 

5 http://wfirst.gsfc.nasa.gov/ 
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best use of WL data in upcoming surveys are still un¬ 
der debate. The problem originates from the fact that 
cosmic shear follows non-Gaussian probabilit y distri bu- 
tion_due to non-linear gravitational growth (Sato et al. 
2009). In order to incorporate the non-linear features 


accurately into WL statistics, cosmological 7V-body sim¬ 
ulations have been extensively used. Previous num erical 
studies ( Hilbert et al.||2009 Sato et ah|2009 2011) have 
already provided important guides for cosmological stud¬ 
ies with WL statistics. There still remain several possible 
factors to be examined. One of the uncertainties in such 
studies is the effects of baryonic physics. Modeling bary- 
onic effects on the WL statistics is difficult because of 
the overall complexities in galaxy formation. Recent nu¬ 
merical simulations and semi-analytic methods success- 


Martizzi et al. 2012 

Schaye et al. 2014; (Jkamoto ef al. 

20141 IMartizzi et al. 

2014j 

Schaller et al. 2014a|b[ Veilis- 


focus on active galactic nuclei (AGN) feedback, which 
quenches star formation in massive halos and may solve 
over cooling problem (i.e. the over production of stars 
in numerical simulations). These simulations also show 
that baryonic physics can change significantly the distri¬ 
bution of both dark matter and baryons within a halo. 
However, since baryonic effects are expected to be weak 
at large scale, most of WL studies so far are based on 
simulations with dark matter component only. For fu¬ 
ture ‘precision cosmology’ with WL, neglecting various 
astrophysical processes may lead to undesirable bias of 


out by several studies (Jing et al. 2006; 

Semboloni et al. 

2011 

Zentner et al. 2013; Yang et al. 2 

1113; MolrammeH 

et al. 

2014J). Hence, it is crucial and timely to study the 


In this paper, we use several statistics to extract the 
non-Gaussian information of WL maps. In addition to 
power spectrum (PS) of weak lensing convergence, we 
conside r peak counts and the Minkowski Functionals 


(MFs). Kratochvil et al. (2010) show that peak counts 
on lensing map can be indeed useful for cosmological pa¬ 
rameter estimation. High-o lensing peaks are likely as- 
sociat ed with massive dark matter halos along a line of 
sight dHaina.ua et al. 2004 Yang et al. 2011 Hamana 


et al. 2012]T MFs are morphological statistics for a given 


multidimensional field. MFs are one of the good mea¬ 
sure of topology in WL and contain the suitabl e cos- 


mological information beyond two-point statistics (Mun- 
slii et aL]|2012t ( Kratoch yil et al. j [2012[ |Shirasa ki et al. 


2012 Shirasaki & Yoshida 201 4 [ ) . Throughout this pa- 


per 7 we use the terms “local statistics” and “non-local 
statistics”. Shear two-point correlation functions and 
the corresponding power spectrum are measured directly 
from shear that is obtained from local measurements, 
i.e., from galaxy ellipticities. Peak counts and MFs are 
non-local statistics, because they can be obtained from 
convergence maps that are derived by integrating shear 
over the region of interest. In order to clarify the bary¬ 
onic effects on the WL statistics, we utilize a large set 
of hydrodynamical ./V-body simulations including vari¬ 
ous baryonic processes. We then study how the baryonic 
effects bias cosmological parameter estimation with WL 
measurement. The rest of this paper is organized as fol¬ 
lows. In Section[2j we summarize the basics of WL statis¬ 


tics of interest and the implementations to estimate the 
statistics for a given WL data. In Section |3j we describe 
our simulation set that includes cosmological simulations 
with baryonic physics. In Section |4j we provide the de¬ 
tails of analysis performed in this paper. In Section [5j 
we show the impact of baryonic physics on WL statis¬ 
tics. We also perform a Fisher analysis to present the 
expected cosmological constraints in upcoming lensing 
surveys. We show the results of \ 2 analysis to quantify 
the baryonic effects on parameter estimation with WL 
statistics. Concluding remarks and discussions are given 
in Section [G] 


2. WEAK LENSING STATISTICS 

We first summarize the basics of gravitational lensing 
by large-scale structure. Weak gravitational lensing ef¬ 
fect is characterized by the image distortion of a source 
object by the following 2D matrix: 

where the observed position of a source object is denoted 
by 0, the true position is f3, k is convergence, and 7 
is shear. In weak held, where gravitational potential is 
small compared with c 2 , each component of Ajj can be 
related to the second derivative of the gravitational po¬ 
tential $ as 


_ d/3 l _ /1 - K - 71 

lj ~ d 0 i ~ V -72 


A-ij — ^ij ^ ij 1 


( 2 ) 


= ^ J q dx7(X,XO^^$KX')0,X']> (3) 


f(x,x') = 


r(x ~ x'Mx') 

r{x) 


(4) 


where r(x) is angular diameter dis tance, an d 37 =_r 0 i 


2001 Munshi et al. 

2008 

). By using the Poisson equa- 

tion and Born approximation (Barteimann &; Schneider 

2001 

Munshi et al. 2008 

, one can express weak lensing 


convergence held as 

, a \ 3 ( H ° 

= O — 


fin 


, d */(X.X) aM <5) 


Born approximation 

yields sufficiently 

accurate two- 

point statistics (e.g., 

Schneider et al. 

1998 

). In the 
inearity of 

present paper, we take into account the non- 


convergence shown in Eq. ([ 3 ]) by performing ray-tracing 
simulations through the matter density held obtained 
from cosmological simulations. 


2.1. Observables 

Here, we summarize three different statistics of weak 
lensing convergence held, power spectrum (PS), peak 
counts, and Minkowski Functionals (MFs). PS has com¬ 
plete cosmological information only if statistical prop¬ 
erties of matter fluctuation follows Gaussian distribu¬ 
tion. However, non-linear structure formation induced 
by gravity inevitably makes the fluctuation deviate sig¬ 
nificantly from Gaussian. In order to extract cosmologi¬ 
cal information, we will also use non-local statistics, i.e. 
peak counts and MFs. 
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2.1.1. Power spectrum 

PS is one of the basic statistics for modern cosmology. 
For a convergence field k. PS is defined by the two point 
correlation in Fourier space: 

<«(£)«* (*)> = (2 (6) 


where the multipole £ is related with angular scale 


through 9 = tt/£. By u sing Limber approximation (Lim¬ 
ber 1954; Kaiser 1992) and Eq. © , one can derive the 


convergence power spectrum as 


p k W 



dx 


w ( x ) 2 

r(x ) 2 



(7) 


where P$(k ) represents the three dimensional matter 
power spectrum, \ a is comoving distance to source galax¬ 
ies and W (x) is the lensing weight function defined as 


w (x) = U—) n 


r{Xs ~ X'Mx) 

r(Xs) 


(1 + Z(X))- (8) 


We follow Sato et al. (2009) to estimate the conver¬ 
gence PS from numerical simulations. We measure the 
binned power spectrum of convergence field by averag¬ 
ing the product of Fourier modes |k(£)| 2 . We use 20 bins 


logarithmically spaced in the range of £ = 100 to 10 5 . 
For parameter estimation performed in section 4, we re¬ 
compute PS using 10 bins logarithmically spaced in the 
range of £ = 100 to 2000. 


2.1.2. Peak count 

Pea ks in convergence maps can be a probe of massive 
halos (Hamana ct al. 2004, 2012) and thus contain cos¬ 
mological information( Yang et al.||2(JII ). 

In practice, peaks on convergence map is defined by a 
local maxima on the “smoothed” map. We do the map¬ 
smoothing because the observed lensing field is signifi¬ 
cantly contaminated by the intrinsic cllipticities of source 
galaxies. The contaminant is called shape noise, which is 
indeed the major contribution to the measured shape of 
source galaxies. In practice, we use a Gaussian window 
function in order to reduce the effect of shape noises on 
WL statistics. The smoothed convergence K.(0\9g) can 
be written as convolution with a filter function of Wq '■ 


by including the statistical properties of shape noise and 
the impact of shape noise on peak position. 

In order to locate peaks on a discretized map obtained 
from numerical simulations, we define the peak as a pixel 
that is higher than eight neighboring pixels. We then 
measure the number density of peaks as a function of 
KL. We exclude the region within 29 g from the bound¬ 
ary of the map in order to avoid the effect of incomplete 
smoothing. In this paper, we divide peaks into two sub¬ 
groups: “medium peaks (MPs)” and “high peaks (HPs)”. 
We define M Ps and HPs by the peak height in a simi¬ 
lar manner to'Yang et al.| (2013). The former is defined 
by the covergence peak with 1.0 < /Cpeak/onoise < 3.0 
whereas the latter corresponds to the peak with 3.0 < 


peak / C 


< 5.0. Here, a, 


5 v noise 


is the rms of shape noise 


on smoothed map given by Eq. (23). It is important to 
note that these peaks are thought to have different phys¬ 
ical origins. MPs are likely caused by the shape noise 
or/an d sever al dark matter halos aligned along a line of 
sight (|Yang et al.|2011 ), whereas HPs are associated with 


indiv idual massive dark matter halos (e.g., Hamana et al. 


20041. Throughout the paper, we set the number of bins 


to be 10 for parameter estimation when measuring peak 
count of each subgroup. 


2.1.3. Minkowski Functionals 

MFs are morphological descriptors for smoothed ran¬ 
dom fields. There are three kinds of MFs for two- 
dimensional maps. Each MFs of Vo, Vp, and Vi repre¬ 
sent the area above the threshold v , the total boundary 
length, the integral of geodesic curvature along the con¬ 
tours. They are given by 


V 0 (v) 


VM 


Vi{v) 



( 11 ) 

( 12 ) 

(13) 


IC(6-9 g ) = j d 2 d> W G (0-0;0 G )/<0), 


(9) 


where 0 G is the smoothing scale and Wq is a gaussian 
filter given by 


Wo(e) = 


( 10 ) 


One can evaluate the smoothed convergence due to an 
isolated massive cluster at a given redshift by assum¬ 
ing the uni versal matter density p rofile of dark matter 


halos (e.g., Navarro et al. 1997b). A simple theoreti¬ 
cal framework to pr edict the number density of the K, 
peaks is presented by Hamana et al. (2004|; Maturi et al. 
(2011). Their calculation yields reasonable results when 
the signal-to-noise ra tio of K, due to mas sive halos is 
larger than ~ 4 (see, Hamana et al. 2004 for details). 


Fan et al. (2010) also consider more detailed calculation 


where K is the geodesic curvature of the contours, da 
and d£ represent the area and length elements, and A 
is the total area. Q u and dQ u are denoted to be excur¬ 
sion sets and boundary sets for the smoothed field IC(x), 
respectively. They are defined by 


Q v = {x\K,(x) > v}, (14) 

dQ v — {x\K(x) = ik}. (15) 


In particular, Vi is equivalent to a kind of genus statistics 
and equal to the number of connected regions above the 
threshold, minus ones below the threshold. Therefore, 
for high thresholds, Vi is essentially equivalent to the 
number of peaks. 

For a two-dimensional Gaussian random field, the ex¬ 
pectation values for MFs can be described by analytic 
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functions as shown in (Tomita 1986): 


Vo{v)= 2 


1 - erf 


v — li 
a Q 


rr , , 1 Ol 

Vi(i') = —7=— exp 
8^2 <r 0 




I/-/Z crj 

F2M = 2(2^^ eXP 


(^-/T 


(16) 

(17) 

(18) 


where // = (/C), cTq = (/C 2 ) — /x 2 , and cr 2 = (|V/C| 2 ). Al¬ 
though MFs can be evaluated p erturbativelv if the non - 
Gaussianity of the field is weak ( Matsubara|2003[ 2010), 
it is difficult to ev aluate MFs of highly non-Gaussian field 
(Petri et al. J20i3). In this paper, we pay a spatial atten¬ 
tion to the non-Gaussian cosmological information ob¬ 
tained from convergence MFs. 

For discretized K, m aps, we employ followin g estima¬ 
tors, as shown in, e.g., Kratochvil et al. (20121, 



V 0 {u) = j J 0 (/C - v)dxdy, 

Vl ^ = IA J 5 ( K '- v ^\I K 'x + dxdy > 

J 6{K. - v) 


Fig. 1.— Top panel : We plot the halo mass functions at z = 0 
measured from our set of simulations. The horizontal axis repre¬ 
sents total mass of halos which include dark matter, gas and star 
particles. The blue, green and red line corresponds to simulation 
results averaged over 10 realizations of DM model, BA model and 

(19) FE model, respectively. The error bars represent the standard de¬ 
viation of the ten realization s for each model . The cyan line shows 
model prediction by|Sheth & Tormen| (|1999|). Bottom panel : The 

(20) fractional difference ol baryomc models trom the fiducial model. 


V 2 (u) = 


2 K X KyK X y - IC 2 x K.yy - 


2nA 


ICl + Kl 


-dxdy, 

( 21 ) 


where 0(x) is the Heaviside step function and 5{x) is 
the Dirac delta function. The subscripts represent dif¬ 
ferentiation with respect to x or y. The first and sec¬ 
ond differentiation are evaluated with finite difference. 
We precompute MFs for 100 equally spaced bins of 
z/ = (v — (/C))/cr n between —10 to 10. For cosmological 
parameter estimation, we recalculate values on equally 
spaced 10 bins in the range —3.0 < v' < 3.0 from 100 
bins. 

3. SIMULATION 
3.1. N-body simulations 

We are interested in non-linear gravitational evolution 
of large-scale structure. In order to follow the evolution 
accurately, we run cosmological iV-body simula tions. We 
use parallelized tree-PM code Gadget-3 (Springel [2005) 
with baryonic processes (discussed below) to follow struc- 
ture formation from an early epoch (z = 99) to present 
(z = 0). The i nitial condi tions are generated by MUSIC 
code (Hahn & Abel 2011), which is based on the sec- 


ond order L agrangian perturbation theory (e.g., Crocce 
et al. 2006). The transfer f unction is generat ed by the 
linear Boltzmann code CAMB (Lewis et al. 2000). The vol¬ 
ume of each simulation box is comoving 240 Mpc /h on 
a side. We adopt the fiducial cosmological parameters 
as follows: matter density H m = 0.279, baryon density 
fib = 0.0463, dark energy density Ha = 0.721, Hub¬ 
ble parameter h = 0.70, spectral index n s = 0.972 and 
amplitude of scalar perturbation A s = 2.41 x 10 -9 at 
the pivot scale k = 0.002 Mpc -1 . These parameters are 
consistent with 9-year Wilkinson Microwave Anisotropy 
Probe (WMAP) result ( Hinshaw et al.||2013j) . 

In this paper, we perform three kinds of cosmologi¬ 
cal simulations; cold dark matter (CDM) simulation (de¬ 


noted as “DM”) and two baryonic simulations. In or¬ 
der to model the degeneracy between cosmological pa¬ 
rameters in WL statistics, we run six CDM only sim¬ 
ulations with one cosmological parameter varied. Cos¬ 
mological parameters for our models are summarized in 
Table [l] The number of particles is set to be 512 3 for 
CDM only simulations and 2 x 512 3 for baryonic sim¬ 
ulations, which consist of both CDM and gas particles. 
For the fiducial cosmological model, the mass of a par¬ 
ticle is found to be m p = 7.97 x 10 9 M Q /h and m c dm = 
6.65 x lO 9 M 0 //i, m gas = 1.32 x 10 9 M 0 //i for baryonic 
simulations. Our baryonic simulations are based on two 
models denoted as BA and FE. BA model contains adi¬ 
abatic gas particles, which exert only adiabatic pressure 
in a Smoothed Particle Hydrodynamics (SPH) manner. 
In FE model, we employ the galaxy formation ‘recipes’ 
of Okamoto et al. (2014). Our FE model corresponds to 
SIS + AGIN model in their paper with some modifications, 
which include star formation, radiative cooling, super¬ 
nova (SN) feedback, stellar wind feedback and an ad hoc 
active galactic nuclei (AGN) feedback. This prescription 
saves computation time and contains less free parameters 
than solving fully the evolution of black holes. In their 
model, the radiative cooling rate is exponentially sup¬ 
pressed when the velocity dispersion within a dark mat¬ 
ter halo exceeds some threshold. Our simulations em¬ 
ploy a simpler method that switches off radiative cooling 
when the local velocity dispersion reaches the threshold. 


Okamoto et al. (2014) report that the sudden change of 
the cooling"function leads to somewhat artificial increase 
of the stellar mass, but we expect this modification does 
not affect the final results significantly. 

In order to check basic statistics of our simulations, 
we measure halo mass function for each model and stel¬ 
lar mass function for FE model. We run a friends-of- 


friends (FoF) halo finder SubFind (Springel et al. 2001 
Dolag et al. 20091 to identify halos in simulations. For 
our simulations with baryons, the procedure is different 
from that for CDM only simulation. First, using only 
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TABLE 1 

Cosmological parameters used for simulations. 


Run 

w 

10 y A s 

iZm 

S2a 

<78 

No. of sim. 

No. of maps 

Explanation 

DM 

-1.0 

2.41 

0.279 

0.721 

0.821 

10 

100 

CDM only fiducial model 

BA 

-1.0 

2.41 

0.279 

0.721 

0.821 

10 

100 

CDM and adiabatic gas 

FE 

-1.0 

2.41 

0.279 

0.721 

0.821 

10 

100 

CDM and baryonic processes 

High H m 

-1.0 

2.41 

0.302 

0.698 

0.872 

10 

100 

la higher model 

Low i 1 j-j-) 

-1.0 

2.41 

0.256 

0.744 

0.767 

10 

100 

la lower model 

High w 

-0.8 

2.41 

0.279 

0.721 

0.766 

10 

100 

higher w model 

Low w 

-1.2 

2.41 

0.279 

0.721 

0.860 

10 

100 

lower w model 

High A s 

-1.0 

2.51 

0.279 

0.721 

0.838 

10 

100 

la higher A s model 

Low A s 

-1.0 

2.31 

0.279 

0.721 

0.804 

10 

100 

lcr lower A s model 


Note. — Two parameters (fi m , 10 9 A S ) are varied by ltr value of WMAP nine-year result (jllinshaw et al. 20lTj|) and we also adjust 
Qa accordingly to keep the universe spatially flat. We also show the normalization of linear matter power spectrum trg, which is a derived 
parameter. 



Fig. 2.— Galaxy stellar mass functi ons at z = 0. The cyan 
points show observational estimates from|Baldry et al.|p012). The 
red line illustrates the result from FE simulation averaged over 10 
realizations and the error bars represent the standard deviation. 


dark matter particles, we find FoF groups which con¬ 
tain more than 31 particles. Next, we find gas particles 
for BA model and gas and star particles for FE model 
linked with a dark matter particle which belongs to a 
FoF group. Finally, we separate each FoF group into a 
central halo and substructures. In FE model, the stel¬ 
lar component of a halo is called ‘galaxy’ and ‘halo’ in¬ 
cludes the baryonic components and dark matter. Fig¬ 
ure [T] and [2] represent halo and stellar mass functions 
measured from simulations a,t z = 0. We com pare the 
simulation resu lt with theoretical prediction of |Sheth~fc 


Tormen| [19991). Figure | 
results are well describe 


[I] clearly shows our simulation 
:a by the theoretical prediction 
over the typical cluster mass scale, while there are ap¬ 
preciable discrepancies between DM and FE simulations 
in the halo mass of less than ~ 10 13 " 5 Mq. This is partly 
induced by SN explosions and/or stellar winds. Energy 
released by SN explosions or stellar winds expel gas parti¬ 
cles from the halo but massive halos can retain the parti¬ 
cles by the deep gravitational potential well. As a result, 
the mass function for our FE model is slightly smaller at 
small masses. We also compare the stellar mass functions 
of our simulation results and the obse rvati on of galaxy 
stellar mass function in Baldry et al. ([2012). The stel¬ 
lar mass distribution in our FE' model is in reasonable 
agreement with the observation at z ~ 0. Because we 
a dopt a simpler model t han the original implementation 
of Okamoto et al. (20141 that includes radiation pressure 
feedback, our FE run produces less small-mass galax¬ 


ies. Also, our simulation results are inaccurate at stellar 
masses less than ~ 10 lo M Q , because such small galaxies 
contain only a few star particles. These small discrep¬ 
ancies are, however, unimportant in the present paper 
because we draw our main conclusions through compar¬ 
isons of multiple sets of simulations with and without 
baryonic components, rather than through comparisons 
of different feedback models. 


3.2. Ray-tracing simulations 

For ray-tracing simulations of gravitational lensing, we 
utilize multiple sim ulation boxes to generate light-cone 


outputs si milarl y to|White fe Hu|( |2000D , |Hamana & Mel-| 
lier (|2001|), and Sato et al.| [2Ul)9|). Details of the config- 
uration arc found in the last reference. 

We place the simulation outputs to fill the past light- 
cone of a hypothetical observer with an angular extent 
5° x 5°, from z = 0 to 1. The angular grid size of our maps 
is set to be 5°/4096 ~ 0.073 arcmin. We randomly rotate 
and shift the simulation boxes in order to avoid the same 
structure appearing multiple times along a line-of-sight. 
In total, we generate 100 independent lensing maps for 
the source redshift of 2 source = 10 We show an example 
of convergence maps obtained from each baryonic model 
(DM, BA, and FE) in Figure [3] Note that each realiza¬ 
tion of three models use the same random seed when we 
generate initial conditions and multiple planes. 

In order to make the mock lensing maps more realistic, 
we add random gaussian noises as shape noise to the 
sim ulated convergence dat a (e.g.,|Kratochvil et al.|[2010|) 
and Shirasaki et al. ( 2Q12|) ): 


(^noise (*^1 5 *^2)Nnoise (*^1 5 *^2 )) 


Ugal^p 


^X2 ,X f 2 5 


( 22 ) 


where 8 Xty is the Kronecker delta symbol, n ga i is the num¬ 
ber density of source galaxies, cr 7 is the rms of shape noise 
and Api x is the solid angle of a pixel. In the following, 
we adopt <r 7 = 0.4, 7T ga i = 30 arcmin -2 . These values are 
expected to be typical for HSC surv ey. 

As described in Section |2.1.2| and |2.1.3| smoothing is 
required to measure non-local statistics of convergence 
from noisy WL data. For a given smoothing scale 9g, 
the rms of the shape noise after gaussian smoothing is 


2 For z S ource = 1, the lensing weight function W(x) has a max¬ 
imum at z ~ 0.5. This means that our lensing statistics are not 
sensitive to the large scale structure above z = 0.7 — 0.8. 
















































































Fig. 3.— We show sample convergence maps of DM, BA and FE simulations without noise and smoothing. Red or blue area corresponds 
to high or low convergence respectively. There are many red points, which are identified as peaks. 


TABLE 2 

The configuration of bins for each observable 


Observable 

Range 

No. of bins 

PS 

100 < l < 2000 

10 

MP 

1.0 ^ A^peak/^noise ^ 3.0 

10 

HP 

3.0 ^ ^--peak/^noise ^ 5.0 

10 

MFs 

-3.0 <{y- (/C»Ar 0 < 3.0 

10 for each 


Note. — We summarize the configurations of our statistical 
analysis. Bins are linearly separated except for PS that are evalu¬ 
ated using logarithmically spaced b ins, 
given by 


van Waerbeke (2000 


= 0.0291 ~~ 


(ZL) 

V0.4/ 


0 G y 


1 arcmin 


ll'gal 


30 arcmin 


-2 


- 1/2 


(23) 


DM models with seven different cosmologies. Thus, we 
obtain the theoretical model of W for a given cosmolog¬ 
ical model as follows: 


Ni(p) = Ni(p°) + £ 

a 


Nito% + ) - Ni(P5r) 


p% + 


(Pa 


Pi), 


(25) 


where a runs 1 to 3 and the superscript 0 means fiducial 
value p° = (0.279,-1.0,2.41). p a± represents a vector 
with one parameter with a higher or lower value; for 
example p 1+ = (0.302,-1.0,2.41) (the other parameter 
values are given in Table (!]). 

The covariance matrix of the data vector TVj on 25 deg 2 
maps is estimated as 


4. ANALYSIS 

We quantify the effects of the baryonic processes on 
weak lensing analyses in terms of errors and bias in cos¬ 
mological parameter estimation. We consider primarily 
a lensing survey with a sky coverage of 1400 deg 2 , i.e., 
the ongoing wide-field survey by Subaru Hyper Suprime- 
Cam (HSC). In the following, we describe in detail the 
calculation of the ensemble averages of three statistics 
and the covariance matrix in order to derive statisti¬ 
cal implications and to estimate cosmological parameters 
p = (f10 9 A S ). 

4.1. Theoretical model and covariance of lensing 
statistics 

First, we calculate the theoretical template and covari¬ 
ance of the WL statistics using our ray-tracing simula¬ 
tions. The configuration of bins are given in Table [2j 
The model template and covariance are based on our 
DM models. Later, we examine the baryonic effects on 
WL cosmological analysis by comparison with DM model 
and two baryonic models. 

We represent the data vector as A/ and denote the 
dimension of the data vector as n. Note that n can be 
larger than ten, when multiple observables are combined. 
We derive theoretical prediction of lensing statistics by 
averaging over R = 100 realizations: 

1 R 

(Ni{p)) ^ Ni{p)= -^Ni(r,p). (24) 

r=l 

In order to consider the cosmological parameter depen¬ 
dence on N{, we employ linear interpolation based on 


1 R 

c ij ( p ) = [ HAi(DP) - Ni{p)][Nji^p) - Nj{p)]. 

r—1 

(26) 

We ignore th e cos mological dependence oiC(p) (see, e.g., 
Eifle r et al.|2009 ) and hence evaluate the covariance ma¬ 
trix by using the fiducial model, i.e. C(p°). For the 
fiducial HSC survey, we simply scale the covariance ma¬ 
trix by survey area, by multiplying the covariance ma¬ 
trix Cij(p°) in Eq. (26) by a factor of 25/1400. When 
calculating the inverse covariance, we include a debi¬ 
asing correction, the so-c alled An derson - Hart lap factor 
a = (R — n — 2 )/{R — 1) ( Hartlap et al. 112007 ), where R 
is the number of realizations and n is the dimension of 
the data vector. 


4.2. Fisher analysis 

Fisher analysis gives a simple forecast for statisti¬ 
cal confidence level of three cosmological parameters 
(fl m , w, 10 9 H S ) with WL statistics. The Fisher matrix 
is given by 

Fi j = \ r & [AiAj+C-'Mij], (27) 


where A = C~ 1 dC/dpi, Mij = 2(dN/dpi)(dN/dpj) 1 . 
The first ter m vanishes w hen the cosmological depen¬ 
dence is weak fEifl er et al.|2009 ). To compute the second 
term, the first derivative is evaluated by the first order 
finite difference, which is given by 


dN 

dpi 


n( P n - jv(pd 


i- 

Pi ~ Pi 


(28) 
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Then, the marginalized error over the other two param¬ 
eters is given by 


= \/( F 1 )i»- 

4.3. Mock tensing surveys 


(29) 


In order to make a realistic forecast in upc omin g HSC 
survey, we employ the bootstrap method as in Yang et al. 
(2013). Since our suite of WL maps consist of one hun- 
dred 25 deg 2 maps, we randomly choose 1400/25 = 56 
realizations from one hundred 25 deg 2 WL maps. Re¬ 
peating this procedure one thousand times, we can get 
one thousand ‘mock’ HSC WL maps. The resulting maps 
are based on the fiducial cosmological model but we also 
have the same set from simulations with the different 
baryonic effects. Lensing statistics of interest (i.e. PS, 
peak counts and MFs) in HSC surveys are evaluated by 
averaging each statistics on a 25 deg 2 map over 56 real¬ 
izations. 

We then utilize the HSC observables derived in this 
way to investigate the impact of baryonic effects on pa¬ 
rameter estimation. We perform y 2 minimization to 
-W.hsc as follows: 


y 2 (r, p, to) = E AN i ( C 1 ) ij AN j > 


(30) 




where ANi = -^ThscM — ^i(p)> r is the index of re¬ 
alizations and 1 < r < 1000, A/ represents the the¬ 
oretical model of A/ for our dark matter only (DM) 
model, and to represents the difference of model of bary¬ 
onic physics, i.e., BA and FE. Suppose we calculate 
the i-th data A/.hsc from a dark matter only conver¬ 
gence map (m = DM), we should then find the result¬ 
ing best-fit points distribute around the fiducial point 
p° = (0.279,-1.0,2.41). However, when one considers 
the case of m = BA or FE, the center of the distri¬ 
bution of best-fit values could be biased in parameter 
space of p if baryonic effects induce discrepancies be¬ 
tween A/,hsc and Ni(p°). That means y 2 values cal¬ 
culated from Eq. (301 do not follow y 2 distribution for 


baryonic models because the ensemble mean and the co- 
variance matrix are computed from the fiducial model, 
not baryonic models. 

There is another way to estimate the bias of parameter 
estimation. The parameter biases of the BA and FE 
models from fiducial param eters can be co mputed in the 
following manner (see, e.g., Huterer et al.j|2006): 


SPa = E F afi E^™ ~ ^(P°)KO* 

hj 


dNf 

'dpp' 


(31) 


where F is the fisher matrix given by Eq. (27) and IV) 


represents the average of A r , over 100 convergence maps 
for to = BA and FE. 

In the following, we quantify the baryonic effects on 
parameter estimation with WL surveys by cons ide ring 
the distribution of y 2 in Eq. ( |30| ) and Sp a in Eq. (31) for 
our BA and FE models. WeTnen compare the olas of 
parame ter estimations with the marginalized error given 
in Eq. (29). 


5.1. Baryonic effect on convergence statistics 

We summarize the main results of the baryonic effects 
on WL statistics in Figure [4j [5] and [6] 

Figure [4] show PS from our maps, where DM, BA and 
FE models are indicated by blue, green and red lines re¬ 
spectively. We also show the theoreticalprediction by the 
cyan line which is calculated by Eq. (17]) with the mod¬ 
eling of m atter power spectrum, essentially an modified 
HaloFit, of Takahashi et al. (2012) We plot the shot noise 
contribution in the bottom panel as the dashed magenta 
line. The error bars represent the standard deviations 
over one hundred maps. The top and bottom panel rep¬ 
resents to results from the maps with and without noises, 
respectively. The lower portion in each panel shows the 
fractional differences of BA and FE models with respect 
to the DM model. 


POWER SPECTRUM 

Compared with DM model, gas pressure suppresses 
small-scale {£ > 10000) structures in BA model. On 
the other hand, in FE model feedback processes are effi¬ 
cient at the scales of 2000 < i < 10000 but cooling pro¬ 
cesses enhance the convergence power at the smaller scale 
i > 20000. These features by c ooling and feedback pro ¬ 
cesses are consistently fo und by |Semboloni et al.| (|2011|) ; 
Mohammed et al. (2014). However, the strength of the 
baryonic enhancement and suppression is different, be¬ 
cause of the details of baryonic physics implementations 
and also partially because of the difference of cosmolog¬ 
ical parameters. 


PEAK COUNT 

Figure[5]shows the peak counts measured from the sim¬ 
ulated convergence maps for the three models. We use 
the same combination of line colors as in Figure [4] The 
error bars represent standard deviations from one hun¬ 
dred maps. The lower portion in each panel shows the 
fractional difference of BA and FE models with respect 
to DM model, and the shaded region indicates the Pois¬ 
son error of the fiducial model, i.e. the square root of 
the number of peaks within a bin, In this figure, we em¬ 
ploy 5 bins per unit S/N ratio v/a no i se , the peak height 
of convergence v divided by cr no j sn . For the a na lytic al 
prediction, we adopt the model o f Hamana et al. (20041 
with Sheth-Tormen (Sheth & Tormen[|115II?Jj) mass func¬ 
tion and Navar ro-Frcnk- White dark matter density pro¬ 
file ( Navarro et al.|1996 1997aI. We assume the relation 
bet ween the hal o mass and concentration parameter as 
in (Hamana et al. 2012), 


c(M, z) = 7.26 


M 


(lO 12 M Q /h 


- 0.086 


(1 + z)- u - 7 \ (32) 


In Figure [5] medium peaks (MPs) correspond to the 
peaks with 0.03 A pea k ^ 0.09, while high peaks (HPs) 
are those with 0.09 ;$ /C pea k ~ 0.15. HPs are typically as¬ 
sociated with massive halos with the mass of ~ 10 14 Mq. 
The number of HPs is less affected by baryonic effects 
because the baryonic processes do not cause very strong 
effect to change the number of such massive halos, as 


5. RESULTS 


seen in Figure 
with masses o: 


MPs often originate from multiple halos 
0 12 — 10 13 Mq aligned in the line of sight 
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Fig. 4.— We plot the convergence power spectrum from one hun¬ 
dred convergence maps. The top and bottom panel show the result 
without and with noise and smoothing. T he solid cyan line shows 
the model prediction using HaloFit of|Takahashi et al.| (|2012|). The 
measured powers from simulations are slightly smaller than the re¬ 
sult of HaloFit due to finite box effect. The dashed magenta line 
in the noisy case shows the power spectrum of the noise. In the 
lower portion of each panel, we plot the fractional differences from 
the fiducial model and the blue filled region corresponds to the ex¬ 
pected statistical error of DM model when one measures the power 
spectrum from one hundred 1400 deg 2 maps. Note that this region 
does not correspond to the error bar, i.e. the standard deviation, 
shown in the upper panel that is measured from our sample of 25 
deg 2 maps. And in the top panel, this expected error of power 
spectrum from noise free maps is so small that it does not appear 
clearly. 




_i_i_i_ £ _U_ 

-0.05 0.00 0.05 0.10 0.15 0.20 


^■"pcak 

Fig. 5. — The number of peaks per unit square degree. The top 
and bottom panel shows the result without and with n oise. The 
solid cyan li ne shows the model prediction of peaks by |Hamana| 
|et al.| (120 04b ■ which assumes that the halo mass function is de¬ 
scribed bypheth &; Tormen| l| 1 999[l a nd the density profile of halo s 
follow Navarro-ktenk-White profile (|Navarro et al.||1996| |1997aJ|. 
The lower portion in each panel shows tne fractional dinerence 
from the fiducial model. The blue filled region represents expected 
Poisson error of DM model from one hundred 1400 deg 2 maps. 


useful to calculate t he e xpec tation values for a Gaussian 
random field (Eqs. (p~6])- (18)) in order to understand the 
overall feature of MFs. In Eqs. (16)-(18), the spectral 
moments <To and cj\ can be calculated from the conver¬ 
gence power spectrum as 


direction. Baryonic physics, such as SN and AGNJeed- 
back, reduce the masses of these halos (see Section 2.1.2). 
Hence the baryonic effects would decrease the number of 
MPs, but the difference is largely made unimportant in 
the noisy convergence maps be cause MPs ar e affected 
significantly by shape noise (Yang et ahpOll). 


MINKOWSKI FUNCTIONALS 

Figure [6] show the MFs computed from our maps of 
DM, BA and FE models. Again we use the same color 
for each model as in the previous figures. The MFs Vo, V\ 
and V 2 are plotted from left to right. The upper (lower) 
two rows represent results from noise-free (noisy) maps. 
Panels in the second and fourth row show the difference 
of BA and FE models from the DM model. The er¬ 
ror bars represent standard deviations from one hundred 
maps. We employ 5 bins per unit 1 / — {v — {1C))/a q. It is 




d H 

(2n 


2 i 2 p p K (£)\w G (u,e G )\ 2 , 


(33) 


where W G represents the Fourier transform of Eq. (101 


which is given by exp (—£ 2 0q/ 4) . Note that W G 

decreases exponentially at £^ 20000 x (1 arcmin/0c). 
^Frorn the result shown in Figure |4j we expect that 00 
of BA and FE models is smaller than that of DM model 
but a 1 would be nearly the same. We confirm this notion 
by direct measurement of a 0 and o\ from our 100 maps. 
It is important to test whether the differences between 
the models are attributed to <7o- In the above analy¬ 
sis, we measured MFs in terms of normalized threshold 
v' because, if the difference is largely owing to varia¬ 
tion of (To 1 the normalization would absorb all or most 
of the differences. The difference in the MFs becomes 
indeed slightly smaller with the normalization, but does 
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Fig. 6. — The three MFs averaged over one hundred convergence maps. Panels in the first and the second rows show the results without 
noise whereas the third and the fourth rows are for those with noise. The error bars show the standard deviation of one hundred maps. 
Note that the range of threshold v is different between the noise-free and noisy cases. The blue filled region corresponds to the expected 
statistical error of DM model from one hundred 1400 deg 2 maps. 


not vanish completely. Therefore, the baryonic processes 
apparently affect the MFs. It is necessary to develop fur¬ 
ther analysis beyond Gaussian description to explain the 
difference of MFs between baryonic and fiducial models. 

5.2. Fisher forecast and bias of cosmological parameters 

We show the results of our parameter estimation in 
Figures [9]to[TT] In each figure, we show two-dimensional 
error contours marginalized over the other parameter. 
The blue, green and red dots indicate the best-fit values 
for DM, BA and FE models, respectively. Red dashed 
lines represent the fiducial values. Note that in these 
three figures, the plotted parameter range is adjusted 
according to the size of the forecast circle for the observ¬ 
able. In order to make comparison between the figures 
easy, we show the lcr confidence regions of six observables 
in the same scale in Figure [8] The biases and marginal¬ 
ized errors of parameters are summarized in Table [3] and 
Figure [7] In Figure [TJ green and red points show the 


biases divided by the marginalized errors for BA and 
FE models. To check the accuracy of our code, we also 
present results from the fiducial simulation. We plot re¬ 
sults obtained from K. maps with shape noise. 

In the following, we summarize the baryonic effects on 
parameter estimation with different convergence statis¬ 
tics. 

POWER SPECTRUM 

The top three panels in Figure [9] show the results of 
parameter estimation for three (DM, BA, FE) models 
using only PS. We find only small bias caused by the 
baryonic physics. That is because, with realistic ground- 
based surveys PS can be measured accurately at low £ 
about 2000 at most, which is the maximum multipole 
used in our estimation. The difference between the fidu¬ 
cial and baryonic models appears at the small scale. The 
PS amplitude and the shape at l ^ 2000 are not signifi¬ 
cantly affected by the baryonic physics. Previous studies 
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bias/error bias/error bias/error 


Fig. 7.— The parameter biases divided by marginalized errors. Green and red points correspond to BA and FE models respectively. The 
biases caused when using only one observable are mostly within la. However, in combined analyses, e.g. PS+MP+HP, the resulting bias 
can exceed the marginalized la error. 


TABLE 3 

The marginalized errors and biases 


Data statistics <5Q m 

5w 

510 a A s 

Fisher forecast marginalized ler error 

PS 

0.0275 

0.276 

1.01 

MP 

0.00853 

0.110 

0.0557 

HP 

0.0129 

0.118 

0.0971 

Vo 

0.00678 

0.0711 

0.0831 

Vi 

0.00767 

0.0911 

0.0931 

v 2 

0.00912 

0.115 

0.0718 

Vo + Vi + V 2 

0.00439 

0.0513 

0.0521 

PS+MP+HP 

0.00478 

0.0520 

0.0371 

PS+MP+Vo 

0.00372 

0.0391 

0.0412 

PS+HP+Vo 

0.00547 

0.0456 

0.0662 

Parameter bias of BA 

PS 

0.00104 

0.000285 

-0.0342 

MP 

-0.00280 

-0.0394 

-0.0187 

HP 

-0.00745 

-0.0648 

0.0278 

Vo 

-0.00385 

-0.0471 

0.0435 

Vi 

-0.00349 

0.0147 

0.101 

v 2 

-0.00518 

-0.0479 

0.0169 

Vo + Vi + Vb 

-0.00271 

0.00352 

0.0477 

PS+MP+HP 

-0.00115 

-0.0145 

-0.0134 

PS+MP+Vb 

0.000493 

0.00803 

0.00428 

PS+HP+Vo 

0.00148 

0.0261 

0.00488 

Parameter bias of FE 

PS 

-0.00171 

-0.0172 

0.0139 

MP 

-0.00591 

-0.0616 

0.00536 

HP 

-0.00668 

-0.0610 

0.0347 

Vo 

-0.000951 

-0.00112 

0.0543 

Vi 

-0.00559 

0.00157 

0.120 

v 2 

-0.00182 

0.0135 

0.0342 

Vo + Vi + V 2 

-0.000412 

0.0347 

0.0481 

PS+MP+HP 

-0.00594 

-0.0587 

0.0156 

PS+MP+Vb 

-0.00288 

-0.00899 

0.0544 

PS+HP+Vo 

0.000589 

0.0136 

0.0241 


Note. — The marginalized errors of three parameters based on 
Fisher forecast of DM model and the biases of BA and FE models 


compared with the DM model. Note that these values are not 
normalized by the fiducial errors. 


( Yang et al.|2013 Mohammed et al.|2014 ) also 
differences due to the choice ol the maximum 


examined 

multipole. 


In our result, the difference at large scale is quite small 
whereas noise dominates at smaller angular scales, where 
PS can not be measured accurately. For this reason, we 
fix the maximum multipole at t = 2000. Then the baryon 
effects on PS is negligible in cosmological parameter es¬ 
timation using PS. 


PEAK COUNT 

The medium (bottom) three panels in Figure [9] show 
the results using only MPs (HPs). As we have discussed 
above, medium height peaks are also dominated by in¬ 
trinsic shape noise but still have cosmological informa¬ 
tion. High peaks have nearly one-to-one relation with 
massive halos, whose masses are less affected by the bary- 



Fig. 8. — We compare the la confidence regions calculated from 
Fisher forecast using six observables separately. The confidence 
region from PS (green line) is clearly the largest. The other non¬ 
local statistics yield substantially smaller confidence regions. 


onic processes. Interestingly, we find that both MPs 
and HPs induce biased parameter estimation of w with 
a level of ~ 0.56cr and 0.52cr, respectively, for FE model. 
The bias may be partly originated from the baryonic ef¬ 
fect on the shape of massive halos. Massive halos tend 
be ro under w hen baryonic processes are included (e.g., 


Kazantzidis et al. 2004). Then the height of a peak asso- 


ciated with a single halo is reduced, and so may be the 
peak counts at K. ~ 0.0 — 0.05. 


MINKOWSKI FUNCTIONALS 


Figure 10 shows the results using each of the three 
MFs, Vo,Vl,V 2 . The center of the dots is shifted in the 
respective parameter space; significant parameter bias 
can be caused by the baryonic effects. Note that the 
absolute shift itself is not very large, but that the bias 
with respect to the error circle is appreciable (see Fig. 
8 for the relative size of the error circles). When one 
constructs the theoretical template of MFs without the 
modeling of baryonic physics, analyses using Vo, Vi and 
V 2 cause the biased estimation of w with a level of ~ 
0.016cr, 0.017cr, and 0.12cr, respectively, for FE model. 
Although it is difficult to explain the origin of this bias 
completely, we expect the following two effects can be 
responsible for the biased parameter estimation: (i) the 
change of variance of K. (i.e. cr 0 ) and (ii) the change 
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Medium peaks 



High peaks 



Fig. 9. — Error contours and best-fit points using PS and peak counts. Blue, green and red dots indicate the best fit values for each 
realization of our DM, BA and FE model suites. Light blue , blue and deep blue contour show lrr, 2<r and 3 a confidence boundary obtained 
in the Fisher analysis of DM model (described in Sec. |4.2| |, respectively. The discrepancy between the distribution of green or red points 
and the blue confidence regions is attributed to “baryonic effects”. Red dashed lines represent fiducial parameter values. 


in halo shape. More detailed modeling with an analytic 
halo approach would be useful to investigate the relation 
between the property of dark matter halos and lensing 
MFs. This is along the line of our ongoing study using a 
large set of cosmological simulations. 

COMBINED ANALYSIS 

Let us consider a combined analysis with multiple ob¬ 
servables in order to tighten the error and possibly mit¬ 
igate the bias due to baryonic effects. 

Our combination of the observables are of the following 


four types: (i) all the MFs, (ii) PS and peak counts, (iii) 
PS, MPs and Vo, and (iv) PS, HPs and Vo- The last two 
are examples of less biased combinations of the statistics. 
We propose the two combinations that are expected to 
cause small net bias of parameters on the basis of the 
result of Eq. (31). The basic idea is to find a combination 
of observables with biases in the opposite directions in 
parameter space. 

Figure[lT]shows the results of the above combined anal¬ 
ysis. All of combined analysis presented here can tighten 
forecast errors, i.e. they can effectively extract more cos- 
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V 0 













Fig. 10. — Error contours and best-fit points using MFs. 

mological information. Some combinations, e.g. PS and 
peak counts, is largely affected by biases by baryonic ef¬ 
fects. This degrades parameter constraints of combined 
analysis. However, if we adopt, e.g. PS+HPs+Vq, such 
combination can mitigate the bias by baryonic effects. 
When all the statistics cause biases with the same sign 
(e.g., for 10 9 A S ), it is safe to combine those with small 
biases. On the other hand, there is still a possibility of 
causing large biases by combined analysis, if one include 
statistic(s) whose bias is very large. In general, cosmo¬ 
logical parameters are degenerated with each other and 
thus it is not trivial to determine the best combination for 
multiple parameters. The above case focusing on a single 


may serve as a useful guide for combined analysis when 
parameter degeneracy is not strong. We further discuss 
with a simple example with details in multi-dimensional 
parameter space of the formalism in Appendix. 

6. CONCLUSION AND DISCUSSION 

We have studied baryonic effects on WL statistics us¬ 
ing a suite of cosmological simulations that incorporate 
galaxy formation processes. Various baryonic processes 
are implemented in our code, such as gas cooling, star 
formation, and stellar feedback. These processes them¬ 
selves are important subjects of research and the model 
uncertainties with free parameters are still controversial, 
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Fig. 11.— Error contours and best-fit points for combined analysis 
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but such simulations can be used to quantify the bary- 
onic effect, or at least to compare the WL statistics with 
those calculated from dark matter only simulations. We 
focus on cosmological parameter estimation using WL 
statistics. To this end, we made realistic mock observa¬ 
tions by performing ray-tracing simulations through the 
non-linear density fields with the size of survey region 
over 1000 square degrees. 

We have studied three statistics, PS, peak counts and 
MFs, calculated directly from mock lensing maps with 
baryons and compare them with the results from our 
‘fiducial’ dark matter only simulations. The PS deviate 
appreciably at l ^ 7000 from the fiducial model due to 
the baryonic proce sses. This feature is also seen in pre¬ 
vious studies, e.g., Semboloni et al. (2011). The shape 
noise dominates, however, over the baryonic effects at the 
small angular scales and thus the baryonic effects are not 
critical in practice in the analysis using PS. Peak counts 
in convergence maps are also affected because the height 
of a peak is sensitive to the mass and the density profile of 
the corresponding halo. Stellar feedback can effectively 
reduce the mass of small galaxy halos, which results in 
decrease of the number of medium height peaks. High 
peaks are less affected by the stellar feedback effects be¬ 
cause the mass distribution in massive, cluster-size lia- 
l os ar e n ot s ignific antly changed by the stellar feedback. 


Yang et al. (2013) examine the influence on peak counts 
of the enhancement of halo concentration parameter, mo¬ 
tivated by radiative cooling. They show that the number 
of HPs is increased accordingly, while the number of MPs 
is mostly unaffected. We do not see the feature in our 
high peak-counts, probably because our simulations in¬ 
clude not only radiative cooling but also supernova and 
stellar feedback that can compensate, at least partially, 
the condensation of baryons. Finally, MFs are morpho¬ 
logical statistics and thus are promising as a probe of 
the baryonic effects that change the internal matter dis¬ 
tribution of halos. The shape noise again substantially 
affects MFs owing to the effectively increased k- variance. 
Also, the smoother distribution of the gaseous compo¬ 
nents causes the overall feature of the MF smoother, es¬ 
pecially in our BA model (see Fig. 6). 

We have considered three kinds of statistics as probes 
of cosmological parameters. When we use only one of 
them, the parameter bias due to the baryonic effects is 
not significant. Most of the biases are found to be within 
ler error for the fiducial 1400 square degree survey. It is 
expected that, with the upcoming Subaru HSC survey, 
cosmological parameters can be determined without be¬ 
ing significantly compromised by baryonic processes as 
long as a single statistic is used. However, because the 
statistical error itself becomes small, roughly in propor¬ 
tion to the increase of square root of the survey area, 
even a very small bias would become critical for future 
surveys with a half or an all sky coverage. The overall 
bias can appear relatively amplified when a combination 
of a few or more statistics are used because the expected 
error of the parameters becomes small and the bias re¬ 
mains the same. For example, when using both PS and 
peak counts, the parameter bias for w is over la. Clearly, 
such bias needs to be well understood before analyzing 
real observational data. It is also desirable to find a com¬ 
bination of statistics that yields both high precision and 
small bias. If we consider only one parameter, a best way 


would be to combine statistics such that their respective 
biases can cancel each other (see Figure 7). Unfortu¬ 
nately, it is generally non-trivial to estimate parameter 
bias in a multi-dimensional space because the degener¬ 
acy between multiple parameters leads to complicated 
dependence on parameters. We suggest to try all possi¬ 
ble combinations and study in detail, as has been done 
in the present paper in the case with only a few statisti¬ 
cal measures. Among the combinations we have tested, 
we find the combination of PS, one of peak counts and 
Vq gives high precision and yet robust results against 
the baryonic effects. Furthermore, all of the quantities 


hashi et al. 

(2012) studied 

approach. 

iamana et al. 


approach is only valid at high S/N ratio ( xi 4) p eaks. 
For low S/N ratio peaks, this relation breaks down. Das| 


& Ostriker (2006) studied the probability distribution 
function of convergence using fitting formula motivated 
by the log-normal distribution. This probability distri¬ 
bution function is di r ectly related with V p . In recent 
studies, | Fed eli (20141; Fedeli et al. (2014); Mohammed 


et al. (2U14p study baryonic effects on PS using a halo 


model. It would be interesting to study the link between 
other two statistics and properties of halos by extending 
the halo model approach. We leave it as a future work. 
The connection between halo properties and WL statis¬ 
tics is also interesting in that we can possibly observe 
halo properties through WL statistics. 

Our results suggest that the non-local statistics have 
a stronger parameter constraining power than PS, even 
when baryonic effects are taken into account. However, 
one needs to be careful in using the non-local statistics. 
Theoretical frameworks to predict PS have been well de¬ 
veloped through both numerical and analytic approach. 
In realistic observations, the shear measurement itself 
is clean and local, especially when compared with the 
other non-local statistics. Overall, PS is the most studied 
statistic and it is intuitively easy to understand based on 



cosmological constraints, there are intrinsic difficulties 
to derive analytical models for MFs and peak counts. 
Hence it is not straightforward to interpret the non-local 
statistics and their dependence on cosmological parame¬ 
ters. In order to rely on the non-local statistics with real 
data, full comprehension of systematic uncertainties in¬ 
duced by observational effects is necessary. Clearly, these 
issues are worth studying further in order to explore fully 
the applicability of the non-local statistics. In this paper, 
we have clarified the theoretical uncertainty of non-local 
statistics due to baryonic processes. Our result is a key 
to the precise application of non-local statistics of cosmic 
shear for cosmological analyses with the unprecedented 
large surveys. 
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Fig. 12.— The parameter biases divided by the square of marginalized errors. Green and red points correspond to BA and FE models 
respectively. 
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APPENDIX 

COMBINED ANALYSIS IN ONE-DIMENSIONAL PARAMETER SPACE 

We present the criterion for the choice of observables for a combined analysis in a simplified case. The parameter 
bias of a combined analysis is given by 




(Al) 


where F and C is the Fisher matrix and the covariance matrix in the case of combined analysis with several observables. 
[)] Assuming that each observable is independent, one can write p a as 




<>Pa ~ {F) a p 


(A2) 

(A3) 


where F s and C s is the Fisher matrix and the covariance matrix for each observable s, e.g. s =PS. 
For one-dimensional parameter space, we can simplify Fisher matrices, 


F s 


-2 


F s 


£n 2 . 


(A4) 


where cr s is the forecast error with the observable s. Hence, the expected bias with a combined analysis is given by 

■sAE-vTe?- (AS) 


When the sum of 6p s /a^ vanishes, the net bias dp also vanishes. Therefore, a combination that makes the sum s mall 
is a good probe with less bias. We show the bias divided by the square of the error for each observable in Figure [12} 
For multi-dimensional parameter space, this proce dure is not generally applicable because the degeneracy between 
different parameters gives additional terms to Eq. (A5). But when the degeneracy is not strong and the bias is not 
large, this condition gives the effective criterion for a less biased combination of observables in terms of cosmological 
parameters. 
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